Outlier Detection & Missing Imputation
Gi-Seop Lee, In-Hyeok Bae
2018-08-24
Data & Customized Function Load
배경 및 목적
최근 많은 분야에서 자동화가 이루어지고 있고, 자연현상의 관측 역시 관측장비의 발달로 자동화가 상당부분 이루어졌다. 생,화학분야는 상대적으로 미흡하지만 조석, 파랑, 수온 등 과 같은 해양 물리자료의 관측은 자동관측 시스템이 대부분이다.
자동화된 관측장비에서 생산되는 시계열 자료는 기존에 사람이 직접 관측하던 것보다 고해상도이기 때문에 기존보다 활용도가 높을 수 있으나, 기술적인 한계로 인한 결측과 이상값, 자연현상의 복잡함을 그대로 기록하는 특성으로 인해 원시자료를 그대로 사용할 때 분석 결과가 왜곡될 위험이 있다.
Acuna and Rodriguez(2004)는 분류모델의 경우, 1%미만: 심각한 영향없음, 1~5%: 관리 가능, 5~15%: 정교한 방법론 필요, 15%이상: 예측 모델의 심각한 성능저하 로 결측치 비율에 따른 결과 왜곡 현상이 있음을 언급하고 있으며, 강필성(2012)은 결측치를 제외하는 방법보다 간단한 선형 대치 기법을 통해서 자료를 보완하는 것이 분류모델의 성능저하가 적다고 언급했다.
이상자료의 경우도 결측과 마찬가지로 자료의 통계량을 계산하거나 모델을 구성할 때 왜곡을 발생시킬 수 있으며(조홍연 등 2016) 자료에 실제 존재하는 극값과의 구분이 어려워 실제 자료분석에서의 이상자료 탐지 및 제거기법의 활용은 미흡하고 경험에 의존하여 제거하기도 한다.
따라서 본 연구에서는 기상청에서 제공하는 해양기상관측부이의 수온자료를 이용하여 결측치의 대치와 이상자료 탐지 및 대치 기법에 대해서 논의한다.
Summary
wt <- kma_buoy_subset$Water_Temperature
## summary data
print(summary_wt <- wt %>% describe)## .
##
## 1 Variables 87672 Observations
## ---------------------------------------------------------------------------
## Water_Temperature
## n missing distinct Info Mean Gmd .05 .10
## 83891 3781 253 1 17 5.438 10.8 11.6
## .25 .50 .75 .90 .95
## 13.1 15.6 20.9 24.0 25.6
##
## lowest : 0.0 0.1 0.2 5.8 5.9, highest: 30.7 30.9 31.3 31.4 31.6
## ---------------------------------------------------------------------------## missing ratio = 4.5%
print(paste0("missing ratio = ", 100 * round(as.numeric(summary_wt$Water_Temperature$counts[2])/as.numeric(summary_wt$Water_Temperature$counts[1]), 4), "%"))## [1] "missing ratio = 4.51%"Raw data plotting
autoplot(kma_buoy_subset$Water_Temperature)
Step 1
Method a: rosner test
가장 먼저 기계적 오류로 인한 이상치로 판단되는 값부터 제거(과정 생략) 일반적으로 특정 값을 선택적으로 제거하거나(ex -999),rosnerTest같은 분포기반 이상치 탐지 기법으로 제거.
rosner test로 제거하는 과정이 계산시간 측면에서 효율적이지 못하여(코드 수행시간이 오래걸린다.) 본 과정에서는 해당 코드를 실행하지 않는다.
Method b: IQR filter
Outlier는 1.5IQR 범위를 벗어나는 mild outlier와 3IQR 범위를 벗어나는 Extreme Outlier로 구분할 수 있음.
여기서는 사전에 정의한 IQR_filter()함수를 이용하여 3IQR(default) 범위를 넘어서는 Extreme Outlier만 탐지하도록 설정.
IQR_filter## function (x, type = "extreme")
## {
## if (type == "extreme") {
## IQR_range <- 3
## }
## else if (type == "mild") {
## IQR_range <- 1.5
## }
## else {
## stop("Error, correct type should be required")
## }
## q <- quantile(x, c(0.25, 0.5, 0.75), na.rm = T)
## lower <- q[1] - IQR_range * IQR(x, na.rm = T)
## upper <- q[3] + IQR_range * IQR(x, na.rm = T)
## return(ifelse(x < lower | x > upper, NA, x))
## }product1 <- kma_buoy_subset$Water_Temperature %>%
apply(2, IQR_filter) %>% as.xts
p1 <- autoplot(product1)
p1 원 자료에 그대로 IQR filter를 적용해도 2009년 근처의 이상치 제거가 되지 않음.
Step 1 Conclusion
Rosner Test는 코드 수행시간이 오래걸리기 때문에 현재 자료는IQR filter()함수 사용이 적합하다고 판단됨.
일괄 처리 문제가 있으므로 Shifting method를 사용.
Step 2
Polynomial interpolation method: NA Interpolation using imputeTS package
Step2-a: Interpolation Method ‘Linear’
Step2-b: Interpolation Method ‘spline’
Step2-c: Interpolation Method ‘Stineman’
Step 2 Conclusion
이상치가 있는 구간에서 결측이 함께 발생할 경우 결측값 추정이 부적절하며,
따라서, 이상치에 강건한
LOWESS(Locally WEighted Polynomial Scatterplot Smoothing)기법을 이용하여 이상치를 탐지 및 제거한다.
또한, Short term missing에 대한 대치 성능은 크게 문제가 없어보이나, 다항식을 이용한 적합 시, 장기결측 구간에서 결측치 대치 성능이 크게 떨어지는 현상이 발생함.
Step 3
Lowess fitting method: Locally Polynomial Smoothing after IQR filter
Interactive plot
LOWESS 적합과정에서 사용된 Bandwidth 는 전체 자료에 대해서는 적절한 Smoothing 으로 판단되지만, 부분별로 확대해서 살펴보게 되면 연별 변화양상보다 짧은 주기의 변동에 대해서는 Over Smoothing 이 된 것을 확인할 수 있다. 따라서 대부분 주기성분을 갖는 해양 관측자료를 Smoothing 할 때에는 가장 짧은 주기의 변동성분을 대표할 수 있는 Bandwidth 설정이 필요하며 특히 장기관측자료를 한번에 처리할 경우 단주기의 정보손실을 고려해야 한다.
아래 코드를 통해 잔차의 특성을 살펴본다.
## residuals1
##
## 1 Variables 87672 Observations
## ---------------------------------------------------------------------------
## Water_Temperature
## n missing distinct Info Mean Gmd .05
## 83891 3781 83891 1 -0.005396 1.292 -1.81137
## .10 .25 .50 .75 .90 .95
## -1.38475 -0.71024 -0.05238 0.66364 1.45776 2.02038
##
## lowest : -24.969102 -19.940401 -19.136007 -19.029932 -19.027409
## highest: 5.961000 6.261774 6.317315 6.418865 6.935617
## ---------------------------------------------------------------------------Step 3 Conclusion
Step 3-a 에서 imputeTS package의 na.interpolation() 함수를 사용하여 Raw data의 결측구간의 보간을 시도함. 짧은 구간의 보간은 비교적 잘 되었으나 장기결측 구간의 보간은 값을 직접 사용하기에는 문제 발생.
따라서 Step 3-b 에서 수행한 Locally Polynomial Smoothing 으로 계산된 값을 이용하여 다음 단계인 Step 4 에서 Step 1의 IQR filter Step 3-a 의 Interpolation과 EMD, SSA 기법을 적용해본다.
Step 4 Iteration
Step 4-a Iterate IQR filter for Smoothed Data to Detect Masked Outlier & Repeat Smoothing
Step 3 에서 계산된 Smoothing Data 에 IQR filter를 다시 적용하며 이는 Masked Outlier 를 탐지해 내기 위함이다. 여기서는 정보손실을 최소화 하기 위하여 Extreme Outlier 탐지 조건인 3*IQR 범위를 적용하였으나, Over Smoothing 된 값으로 인해 Raw Data 의 정상치로 추정되는 값이 이상자료로 판단되어 제거되는 것을 방지하기 위해 더 완화된 3*IQR 이상의 조건을 고려해 볼 수 있다. 이는 자료를 Smoothing 할 때 발생하는 과적합 문제를 해결하기 위해 임의로 Bandwidth를 넓게 조정한 것에 대한 조치이다.
Raw Data 와 Smoothing Data의 잔차를 이용하여 이상자료를 탐지, 제거하고 난 후의 plot.
연별 변화패턴에 Masking 되어 있던 이상자료 중 극단적인 값들이 제거된 것을 볼 수 있으며 이로 인해 다시 계산한 Smoothing line 이 달라짐.
number of outliers removed by ‘residual IQR filter’
## [1] 6740## [1] 254## [1] 281.5 IQR, 3 IQR, 5 IQR일 경우, 각각 6740, 254, 28개의 이상치가 탐지됨.
결측 구간을 제외한 전체 자료 83891개 중 0.3% 미만의 극히 일부자료이기 때문에 본 분석에서는 3 IQR 조건과 5 IQR 조건의 차이는 미미하다고 할 수 있고 따라서 3 IQR조건을 적용하여 이론적인 근거(Extream Outlier)를 확보.
단, 1.5 IQR범위에서 필터링 할 경우, 전체 자료의 약 8%가 이상자료로 감지되기 때문에 고도화된 대치(Imputation)방법 강구 필요.
Step 5 Missing Imputation
Step 1~4 과정을 통하여 Outlier 제거와 Missing Data 의 대치는 동시에 수행되어야 한다는 것을 알 수 있음.
따라서 Step 5 에서는 Step 4 까지의 Outlier 제거과정을 거치며 얻은 Preprocessed Data 를 이용하여 그럴듯 하게 Missing Data Imputation 처리. 각 기법의 적용은 imputeTS package 의 na.kalman(), na.locf(), na.ma(), na.seadec() 함수들과 그 인자를 조정해가며 테스트해본다.
na.seadec() 함수는 자료의 seasonality를 계산하지 못하여 사용불가.